Εργαστηριακή Άσκηση 1#
Σκοπός της πρώτης σειράς ασκήσεων είναι, αφ’ ενός η εξοικείωση με το προγραμματιστικό περιβάλλον της Python, αφ’ ετέρου, η εισαγωγή στους τρόπους παράστασης και επεξεργασίας τηλεπικοινωνιακών σημάτων στη συγκεκριμένη γλώσσα προγραμματισμού.
Μέρος 1: Εξοικείωση με το προγραμματιστικό περιβάλλον της Python#
Λίγα λόγια για τον interactive python interpreter [ipython] και το jupyter notebook
Interactive python interpreter (IPython) notebook#
Το IPython (Interactive Python) (https://ipython.org/) είναι ένα διαδραστικό κέλυφος γραμμής εντολών σχεδιασμένο για τη γλώσσα προγραμματισμού Python.Υπερβαίνει τον παραδοσιακό διερμηνέα Python και προσφέρει βελτιωμένες δυνατότητες (REPL):
1)Ανάγνωσης (Read)
2)Aξιολόγησης (Evaluate)
3)Eκτύπωσης (Print)
4)Επαναλάψης (Loop)
Jupyter notebook#
To Project Jupyter είναι ένας μη κερδοσκοπικός οργανισμός με αποστολή τη συγγραφή ανοικτού λογισμικού για διαδραστικές εφαρμογές. Ξεκίνησε από ipython, ωστόσο σήμερα προσφέρει προγράμματα σε πολλές γλώσσες προγραμματισμού.
Το jupyter notebook (https://jupyter.org/) είναι μια πλατφόρμα web ανοικτού λογισμικού για την ανάπτυξη διαδραστικών εφαρμογών, κυρίως για επεξεργασία (επιστημονικών) δεδομένων και μηχανικής μάθησης.
Εξάσκηση#
from scipy import signal
# Ανατρέξτε στην τεκμηρίωση της βιβλιοθήκης scipy.signal
# https://docs.scipy.org/doc/scipy/reference/signal.html
from scipy.fft import fft, fftfreq
from dash import Dash, dcc, html, Input, Output
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
from numpy import random
import pandas
import plotly
import plotly.express as px
import plotly.graph_objs as go
# Δείτε την έκδοση της αριθμητικής βιβλιοθήκης numpy
import dash_bootstrap_components as dbc
np.__version__
---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Cell In[1], line 5
2 # Ανατρέξτε στην τεκμηρίωση της βιβλιοθήκης scipy.signal
3 # https://docs.scipy.org/doc/scipy/reference/signal.html
4 from scipy.fft import fft, fftfreq
----> 5 from dash import Dash, dcc, html, Input, Output
6 import matplotlib
7 import matplotlib.pyplot as plt
ModuleNotFoundError: No module named 'dash'
import warnings
warnings.filterwarnings('ignore')
Μην ξεχνάτε ότι η IPython μας δίνει τη δυνατότητα να ‘εξερευνήσουμε’ το περιεχόμενο ενός package, χρησιμοποιώντας τη δυνατότητα του tab-completion, ή τη χρήση του ? για help/documentation: Π.χ., για να δούμε όλα τα περιεχόμενα του signal namespace δίνουμε:
In [3]: signal?
και για να καλέσουμε την ενσωμετωμένη τεκμηρίωση της numpy, δίνουμε:
In [4]: np?
Περισσότερες πληροφορίες μπορείτε να πάρετε από το http://www.numpy.org.
# Δημιουργήστε ένα βαθμωτό (μονοδιάστατο) μέγεθος
s=2
print('s=',s)
s= 2
# Δημιουργείστε ένα διάνυσμα πραγματικών τιμών:
# Στο MATLAB: v = [1,5,9] ή v = [1 5 9]
v=np.array([1,5,9])
print('v=',v)
v= [1 5 9]
# Δημιουργείστε έναν πίνακα πραγματικών τιμών:
# Στο MATLAB: a = a=[[1,2,3];[4,5,6];[7,8,9]] ή a=[1,2,3;4,5,6;7,8,9]
a=np.array([[1,2,3],[4,5,6],[7,8,9]])
print('a=',a)
a= [[1 2 3]
[4 5 6]
[7 8 9]]
# Αθροίστε
a+5
array([[ 6, 7, 8],
[ 9, 10, 11],
[12, 13, 14]])
#Πολλαπλασιάστε
b=s*v*2
print('b=',b)
b= [ 4 20 36]
# Πολλαπλασιάστε στοιχείο-προς-στοιχείο (elementwise)
# MATLAB: v.*b
np.multiply(v,b)
array([ 4, 100, 324])
#Ελέγξτε το μήκος ενός διανύσματος
# MATLAB: length(v)
len(v)
3
# Ελέγξτε το μέγεθος ενός πίνακα
# MATLAB: size(a)
a.shape # για array: np.array(a.shape)
(3, 3)
# Προσπελάστε συγκεκριμένα στοιχεία ενός πίνακα
# Η δεικτοδότηση αρχίζει από το 0.
# MATLAB: a(1,2)
# --- ΠΡΟΣΟΧΗ, στο MATLAB η δεικτοδ΄ότηση αρχ΄ίζει από το 1!
a[0,1]
2
# Προσπελάστε συγκεκριμένα στοιχεία ενός πίνακα (συνέχεια)
# Αρνητικές τιμές μετρούν από το τέλος, π.χ. το -1
# αναφέρεται στο τελευταίο στοιχείο)
a[1,-1]
6
# Προσπελάστε συγκεκριμένο τμήμα ενός διανύσματος
# MATLAB: v(1:9)
v[1:3]
# ΠΡΟΣΟΧΗ: τα στοιχεία [2ο,3ο] δίνονται ως 1:3 και όχι ως 1:2
# Δοκιμάστε το v[1:2]...
array([5, 9])
# Προσπελάστε συγκεκριμένα τμήματα ενός πίνακα
a[0:2,:]
# Ομοίως: οι γραμμές 1 & 2 δίνονται ως 0:2 και όχι ως 0:1
array([[1, 2, 3],
[4, 5, 6]])
# Δημιουργήστε ένα διάνυσμα με στοιχεία από το 0 έως το 0.5 και βήμα 0.1
# MATLAB: t=(0:0.1:0.4)
t=np.arange(0,0.5,0.1)
print('t=',t)
t= [0. 0.1 0.2 0.3 0.4]
Μέρος 2: Δειγματοληψία - Ψηφιοποίηση#
Τα πρωτογενή σήματα είναι κυρίως αναλογικά (συνεχούς χρόνου). Για να τα παραστήσουμε και επεξεργαστούμε στον υπολογιστή μας (ή άλλη ψηφιακή μηχανή) θα πρέπει πρώτα να τα ψηφιοποιήσουμε. Υποθέστε ένα σήμα συνεχούς χρόνου \(x(t)\) με μετασχηματισμό Fourier (Continuous Time Fourier Transform – CTFT): $\(X(f)=\int_{-\infty}^{\infty} x(t)e^{-j2\pi ft} dt\)$
app1 = Dash(__name__, external_stylesheets=[dbc.themes.BOOTSTRAP])
# Παράμετροι δειγματοληψίας
fs = 1000 # Sampling frequency in Hz
T = 1 / fs # Sampling period
L = 1000 # Number of samples
# Create time vector
t = np.arange(0, L) * T
# Layout of the Dash app
app1.layout = html.Div(
children=[
html.H2('Interactive Signal and Fourier Transformation'),
dcc.RangeSlider(
id='frequency-slider',
min=1,
max=51,
step=1,
value=[5, 10],
marks={i: f'{i} Hz' for i in range(1, 51, 5)}
),
dbc.Row(
children=[
dbc.Col(
dcc.Graph(id='original-signal'),
width=6
),
dbc.Col(
dcc.Graph(id='fourier-transform'),
width=6
)
]
)
],
style={'backgroundColor': 'white', 'padding': '20px'}
)
# Callback to update graphs
@app1.callback(
[Output('original-signal', 'figure'),
Output('fourier-transform', 'figure')],
[Input('frequency-slider', 'value')]
)
def update_graph(selected_frequencies):
# Generate the signal with the selected frequencies
signal = np.zeros(len(t))
for f in range(selected_frequencies[0], selected_frequencies[1] + 1):
signal += np.sin(2 * np.pi * f * t)
# Calculate the FFT of the original signal
original_signal_fft = fft(signal)
# Calculate frequencies for the FFT of the original signal
original_freqs = fftfreq(L, T)[:L // 2]
# Calculate magnitude of Fourier coefficients (amplitude) for the original signal
original_magnitude = np.abs(original_signal_fft)[:L // 2]
# Original signal graph
original_signal_figure = {
'data': [go.Scatter(x=t, y=signal, mode='lines', line=dict(color='#00CC96', width=2))],
'layout': go.Layout(
title='Original Signal',
xaxis={'title': 'Time (s)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Fourier Transformation graph
fourier_transform_figure = {
'data': [go.Scatter(x=original_freqs, y=original_magnitude, mode='lines', line=dict(color='#1F77B4', width=2))],
'layout': go.Layout(
title='Fourier Transformation of Original Signal',
xaxis={'title': 'Frequency (Hz)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Return this additional figure in the callback
return original_signal_figure, fourier_transform_figure
if __name__ == '__main__':
app1.run_server(debug=True, port=8050)
Λαμβάνοντας δείγματα του \(x(t)\) με ρυθμό \(f_s=1/T_s\) παράγεται σήμα διακριτού χρόνου \(x(nT_s)\). Μαθηματικά το αναπαριστάνουμε ως σειρά συναρτήσεων δέλτα $\(x_\delta (t)=\sum_{n=-\infty}^{\infty}x(nT_s)\delta(t-nT_s)=x(t)\sum_{n=-\infty}^{\infty}\delta(t-nT_s)\)\( με μετασχηματισμό Fourier \)\(X_\delta (f)=\sum_{n=-\infty}^{\infty}x(nT_s)e^{-j2\pi fnT_s}=X(f)*1/T_s\sum_{n=-\infty}^{\infty}\delta(f-k/T_s)=1/T_s\sum_{n=-\infty}^{\infty}X(f-k/T_s)\)$ που είναι περιοδική συνάρτηση.
from scipy.fftpack import fftshift, ifftshift
app2 = Dash(__name__, external_stylesheets=[dbc.themes.BOOTSTRAP])
# Παράμετροι δειγματοληψίας
T = 1 / fs # Sampling period
L = 1000 # Number of samples
# Create time vector
t = np.arange(0, L) * T
# Layout of the Dash app
app2.layout = html.Div(
children=[
html.H2('Interactive Signal and Fourier Transformation'),
dcc.Slider(
id='samples-slider',
min=100,
max=1000,
step=25,
value=500,
marks={i: f'{i} samples' for i in range(100, 1000, 100)},
included=False, # To remove the color fill
updatemode='drag', # Update slider continuously while dragging
tooltip={'placement': 'bottom'} # Show tooltip at the bottom
),
dcc.RangeSlider(
id='frequency-slider',
min=1,
max=51,
step=1,
value=[5, 10],
marks={i: f'{i} Hz' for i in range(1, 51, 5)}
),
dbc.Row(
children=[
dbc.Col(
dcc.Graph(id='sampled-signal'),
width=6
),
dbc.Col(
dcc.Graph(id='sampled-fourier-signal'),
width=6
)
]
)
],
style={'backgroundColor': 'white', 'padding': '20px'}
)
# Callback to update graphs
@app2.callback(
[Output('sampled-signal', 'figure'),
Output('sampled-fourier-signal', 'figure')],
[Input('frequency-slider', 'value'),
Input('samples-slider', 'value')]
)
def update_graph(selected_frequencies,selected_samples):
# Generate the signal with the selected frequencies
signal = np.zeros(len(t))
for f in range(selected_frequencies[0], selected_frequencies[1] + 1):
signal += np.sin(2 * np.pi * f * t)
# Correct downsampling approach
sample_rate = selected_samples # Directly using selected_samples as an integer
downsampling_factor = int(fs / sample_rate) # Downsampling factor
# Select every nth sample from the original time vector to match the downsampled signal
sample_points = t[::downsampling_factor]
# Downsample the signal by selecting every nth sample
downsampled_signal = signal[::downsampling_factor]
# Recalculate L for the downsampled signal if needed
L_downsampled = len(downsampled_signal)
# Calculate the DFT of the downsampled signal
sampled_signal_fft = fft(downsampled_signal)
# Calculate frequencies for the FFT of the downsampled signal
# Note: The new sampling period is the inverse of the new sampling rate
sampled_freqs = fftfreq(L_downsampled, 1/sample_rate)[:L_downsampled // 2]
# Calculate magnitude of Fourier coefficients (amplitude) for the downsampled signal
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
n = len(downsampled_signal)
T = 1 / sample_rate # Recalculate the sampling period with the selected sample rate
# After calculating the magnitude of Fourier coefficients
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
# Create a new frequency vector that includes the replicated frequencies
# First, calculate the original frequency bins for the positive frequencies
original_sampled_freqs = fftfreq(L_downsampled, T)[:L_downsampled // 2]
shifted_freqs_negative = original_sampled_freqs - 1/T # Shift for -1/Ts
shifted_freqs_positive = original_sampled_freqs + 1/T # Shift for +1/Ts
# Concatenate the original and shifted (replicated) frequencies and magnitudes
# This includes the original frequencies, and the replications at -1/Ts and +1/Ts
replicated_freqs = np.concatenate([shifted_freqs_negative, original_sampled_freqs, shifted_freqs_positive])
replicated_magnitude = np.concatenate([sampled_magnitude, sampled_magnitude, sampled_magnitude])
# Sort the replicated frequencies and magnitudes in ascending order for proper plotting
sorted_indices = np.argsort(replicated_freqs)
sorted_replicated_freqs = replicated_freqs[sorted_indices]
sorted_replicated_magnitude = replicated_magnitude[sorted_indices]
# Sampled signal graph
sampled_signal_figure = {
'data': [go.Scatter(x=sample_points, y=downsampled_signal, mode='markers', marker=dict(color='#00CC96', size=7))],
'layout': go.Layout(
title='Sampled Signal',
xaxis={'title': 'Time (s)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Fourier Transformation graph for the sampled signal
sampled_fourier_transform_figure = {
'data': [go.Scatter(x=sorted_replicated_freqs, y=sorted_replicated_magnitude, mode='lines', line=dict(color='#1F77B4', width=2))],
'layout': go.Layout(
title='Fourier Transformation of Sampled Signal',
xaxis={'title': 'Frequency (Hz)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Return this additional figure in the callback
return sampled_signal_figure, sampled_fourier_transform_figure
if __name__ == '__main__':
app2.run_server(debug=True, port=8051)
Για βαθυπερατά σήματα \(x(t)\) εύρους ζώνης W, με την υπόθεση ότι ο ρυθμός δειγματοληψίας \(fs ≥ 2W\), ισχύει ότι \(X(f) = T_s X_\delta(f)\), \(0 ≤ f ≤ W\), δηλαδή, το σήμα \(X(f)\) προκύπτει μετά από διάβαση του δειγματοληπτημένου \(x_\delta(t)\) μέσω ιδανικού βαθυπερατού φίλτρου κέρδους \(T_s\). Από το προηγούμενο σχήμα γίνεται φανερό ότι εάν η δειγματοληψία γίνει με συχνότητα μικρότερη του διπλασίου της ανώτερης συχνότητας \(W\) του σήματος (υποδειγμάτιση – undersampling), τότε εμφανίζονται στην περιοχή συχνοτήτων του σήματος «είδωλα» φάσματος από ανώτερες συχνότητες που δεν επιτρέπουν την ακριβή αποκατάσταση του αρχικού σήματος συνεχούς χρόνου. Το φαινόμενο αυτό ονομάζεται αναδίπλωση ή επικάλυψη (aliasing), το δε σφάλμα κατά την αποκατάσταση του αρχικού σήματος αποκαλείται σφάλμα αναδίπλωσης (aliasing error). Η δειγματοληψία στο πεδίο του χρόνου αποτελεί τη βάση για τον ορισμό του μετασχηματισμού Fourier διακριτού χρόνου (Discrete Time Fourier Transform – DTFT). Για μια σειρά διακριτών αριθμών \(x[n]\), ο μετασχηματισμός Fourier διακριτού χρόνου ορίζεται ως:
O DTFT είναι περιοδική συνάρτηση με περίοδο \(1\), επομένως, αρκεί ο υπολογισμός του στο διάστημα συχνοτήτων \([0,1]\) ή ισοδύναμα \([-½,½]\). Να σημειωθεί ότι ο DTFT, παρότι προκύπτει από μια σειρά διακριτών αριθμών \(x[n]\), είναι συνεχής συνάρτηση της μεταβλητής \(\phi\) όπως παραστατικά φαίνεται στο επόμενο σχήμα.
from scipy.fftpack import fftshift, ifftshift
app3 = Dash(__name__, external_stylesheets=[dbc.themes.BOOTSTRAP])
# Παράμετροι δειγματοληψίας
T = 1 / fs # Sampling period
L = 1000 # Number of samples
# Create time vector
t = np.arange(0, L) * T
# Layout of the Dash app
app3.layout = html.Div(
children=[
html.H2('Interactive Signal and Fourier Transformation'),
dcc.Slider(
id='samples-slider',
min=100,
max=1000,
step=25,
value=500,
marks={i: f'{i} samples' for i in range(100, 1000, 100)},
included=False, # To remove the color fill
updatemode='drag', # Update slider continuously while dragging
tooltip={'placement': 'bottom'} # Show tooltip at the bottom
),
dcc.RangeSlider(
id='frequency-slider',
min=1,
max=51,
step=1,
value=[5, 10],
marks={i: f'{i} Hz' for i in range(1, 51, 5)}
),
dbc.Row(
children=[
dbc.Col(
dcc.Graph(id='sampled-signal'),
width=6
),
dbc.Col(
dcc.Graph(id='sampled-fourier-signal'),
width=6
)
]
)
],
style={'backgroundColor': 'white', 'padding': '20px'}
)
# Callback to update graphs
@app3.callback(
[Output('sampled-signal', 'figure'),
Output('sampled-fourier-signal', 'figure')],
[Input('frequency-slider', 'value'),
Input('samples-slider', 'value')]
)
def update_graph(selected_frequencies,selected_samples):
# Generate the signal with the selected frequencies
signal = np.zeros(len(t))
for f in range(selected_frequencies[0], selected_frequencies[1] + 1):
signal += np.sin(2 * np.pi * f * t)
# Correct downsampling approach
sample_rate = selected_samples # Directly using selected_samples as an integer
downsampling_factor = int(fs / sample_rate) # Downsampling factor
# Downsample the signal by selecting every nth sample
downsampled_signal = signal[::downsampling_factor]
# Recalculate L for the downsampled signal if needed
L_downsampled = len(downsampled_signal)
# Calculate the DFT of the downsampled signal
sampled_signal_fft = fft(downsampled_signal)
# Calculate magnitude of Fourier coefficients (amplitude) for the downsampled signal
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
n = len(downsampled_signal)
T = 1 / sample_rate # Recalculate the sampling period with the selected sample rate
# After calculating the magnitude of Fourier coefficients
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
# Create a new frequency vector that includes the replicated frequencies
# First, calculate the original frequency bins for the positive frequencies
original_sampled_freqs = fftfreq(L_downsampled, T)[:L_downsampled // 2]
shifted_freqs_negative = original_sampled_freqs - 1/T # Shift for -1/Ts
shifted_freqs_positive = original_sampled_freqs + 1/T # Shift for +1/Ts
# Concatenate the original and shifted (replicated) frequencies and magnitudes
# This includes the original frequencies, and the replications at -1/Ts and +1/Ts
replicated_freqs = np.concatenate([shifted_freqs_negative, original_sampled_freqs, shifted_freqs_positive])
replicated_magnitude = np.concatenate([sampled_magnitude, sampled_magnitude, sampled_magnitude])
# Sort the replicated frequencies and magnitudes in ascending order for proper plotting
sorted_indices = np.argsort(replicated_freqs)
sorted_replicated_freqs = replicated_freqs[sorted_indices]
sorted_replicated_magnitude = replicated_magnitude[sorted_indices]
# Generate an array of sample indices for the downsampled signal
sample_indices = np.arange(len(downsampled_signal))
# Update the Sampled signal graph to plot against sample indices
sampled_signal_figure = {
'data': [go.Scatter(x=sample_indices, y=downsampled_signal, mode='markers', marker=dict(color='#00CC96', size=7))],
'layout': go.Layout(
title='Sampled Signal',
xaxis={'title': 'Sample Number', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Normalize the frequency values
normalized_freqs = original_sampled_freqs / (1/T)
# Since you previously concatenated and sorted for replication, ensure to apply normalization there as well
normalized_replicated_freqs = np.concatenate([
(shifted_freqs_negative / (1/T)), # Normalize -1/Ts shifted frequencies
normalized_freqs, # Already normalized original frequencies
(shifted_freqs_positive / (1/T)) # Normalize +1/Ts shifted frequencies
])
# Sort the normalized and replicated frequencies for proper plotting
sorted_indices = np.argsort(normalized_replicated_freqs)
sorted_normalized_replicated_freqs = normalized_replicated_freqs[sorted_indices]
sorted_replicated_magnitude = replicated_magnitude[sorted_indices]
# Update the Fourier Transformation graph with normalized frequencies
sampled_fourier_transform_figure = {
'data': [go.Scatter(x=sorted_normalized_replicated_freqs, y=sorted_replicated_magnitude, mode='lines', line=dict(color='#1F77B4', width=2))],
'layout': go.Layout(
title='Normalized Fourier Transformation of Sampled Signal',
xaxis={'title': 'Normalized Frequency (f/fs)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Return this additional figure in the callback
return sampled_signal_figure, sampled_fourier_transform_figure
if __name__ == '__main__':
app3.run_server(debug=True, port=8051)
Με τη σειρά των διακριτών αριθμών να προκύπτει ως αποτέλεσμα δειγματοληψίας, \(x[n]=x(nT_s)\), ο DTFT και ο μετασχηματισμός Fourier \(X_\delta(f)\) του δειγματοληπτημένου σήματος συνδέονται μέσω της αντιστοιχίας \(\phi ↔ f/f_s\). Η συνήθης πρακτική είναι να παριστάνουμε τον λόγο \(f/f_s\) ως κανονικοποιημένη συχνότητα \(\phi\) (\(f_D\), στις σημειώσεις σας) και οι πραγματικές συχνότητες να προκύπτουν ως πολλαπλάσιά της (συνήθως κλασματικά). Για τη σύνδεση του DTFT με τον μετασχηματισμό Fourier \(X(f)\) του σήματος πρέπει επιπλέον να γίνει αναγωγή στην περίοδο δειγματοληψίας με πολλαπλασιασμό επί \(T_s\) (ή διαίρεση με \(f_s\)). Κατ΄ αναλογία με τη δειγματοληψία σημάτων στο χρόνο μπορούμε να κάνουμε δειγματοληψία στο πεδίο της συχνότητας λαμβάνοντας διακριτές τιμές \(X(kf_o)\) του μετασχηματισμού Fourier που αντιστοιχούν σε ανάλυση συχνότητας \(f_o=1/T_o\). Αυτό ισοδυναμεί με περιοδική επανάληψη του σήματος συνεχούς χρόνου \(x(t)\) κάθε \(Τ_ο\), αφού το περιοδικό σήμα
έχει μετασχηματισμό Fourier
Επομένως, \(X[k] = X(kf_o)/Τ_o\) είναι οι συντελεστές του αναπτύγματος σε σειρά Fourier.του περιοδικού σήματος \(x_p(t)\). Προφανώς, για σήματα \(x(t)\) πεπερασμένης διάρκειας, όπου \(x(t)=0\) για \(|t| ≥ T\), με την υπόθεση ότι η περίοδος \(T_o ≥ 2T\), ισχύει ότι \(x(t) = x_p(t)\) για \(|t| ≤ T\). Στην πράξη, τα σήματα έχουν πολύ μεγάλη διάρκεια για να μπορέσουμε να τα αναλύσουμε στην ολότητά τους. Έτσι εφαρμόζουμε ένα ορθογωνικό χρονικό παράθυρο, ώστε να διατηρήσουμε μόνο το πιο σημαντικό τους μέρος για το διάστημα παρατήρησης και \(x(t)= 0\), αλλού. Κατά τον υπολογισμό του DTFT \(X_d(\phi)\) ενός τέτοιου ακρωτηριασμένου σήματος, αντί του απείρου αθροίσματος, περιοριζόμαστε σε μια πεπερασμένου μήκους \(L\) σειρά αριθμών \(x[n]\), οπότε
H δειγματοληψία του \(X_d(\phi)\) στο πεδίο συχνότητας σε \(Ν\) ισαπέχουσες κανονικοποιημένες συχνότητες \(0\), \(1/Ν\), \(2/Ν\), \(…\), \((Ν-1)/Ν\), δίνει
όπου, εάν \(N≥L\), θέτουμε \(x[n]=0\) για \(n≥L\). Η τελευταία σχέση αναγνωρίζεται ως ο διακριτός μετασχηματισμός Fourier (Discrete Fourier Transform – DFT), ο οποίος για μια πεπερασμένη σειρά \(xn\), \(n=0\), \(1\), \(…\), \(N-1\), ορίζεται ως:
και ο αντίστροφός του είναι
Η \(X_d(\phi)\) ως DTFT είναι περιοδική συνάρτηση και εάν η αρχική σειρά xn ήταν περιοδική (και δεν εφαρμόζαμε το παράθυρο), τότε η \(X_d(\phi)\) θα ήταν μηδέν παντού εκτός των σημείων της δειγματοληψίας \(k/Ν\). Δηλαδή, εάν θεωρήσουμε μια πεπερασμένου μήκους σειρά αριθμών που επαναλαμβάνεται περιοδικά, o διακριτού χρόνου μετασχηματισμός Fourier της (DTFT) είναι και αυτός περιοδικός και διακριτός. Επιπλέον, ο DFT και ο αντίστροφός του IDFT, εάν δεν περιορίζαμε τους δείκτες \(n\) και \(k\) μεταξύ \(0\) και \(N-1\), θα ήταν περιοδικές συναρτήσεις. Άρα η πεπερασμένη σειρά xn μπορεί να θεωρηθεί ως ένα περιοδικό σήμα διακριτού χρόνου ιδωμένο μόνο κατά τη διάρκεια μιας περιόδου και ο DFT, η σειρά \(X_k\), ως τα δείγματα με ανάλυση \(1/Ν\) του DTFT \(X_d(\phi)\) στο πεδίο κανονικοποιημένων συχνοτήτων \([0,1]\), όπως φαίνεται στο επόμενο σχήμα.
from scipy.fftpack import fftshift, ifftshift
app4 = Dash(__name__, external_stylesheets=[dbc.themes.BOOTSTRAP])
# Παράμετροι δειγματοληψίας
T = 1 / fs # Sampling period
L = 1000 # Number of samples
# Create time vector
t = np.arange(0, L) * T
# Layout of the Dash app
app4.layout = html.Div(
children=[
html.H2('Interactive Signal and Fourier Transformation'),
dcc.Slider(
id='samples-slider',
min=100,
max=1000,
step=25,
value=500,
marks={i: f'{i} samples' for i in range(100, 1000, 100)},
included=False, # To remove the color fill
updatemode='drag', # Update slider continuously while dragging
tooltip={'placement': 'bottom'} # Show tooltip at the bottom
),
dcc.RangeSlider(
id='frequency-slider',
min=1,
max=51,
step=1,
value=[5, 40],
marks={i: f'{i} Hz' for i in range(1, 51, 5)}
),
dcc.Slider(
id='N-slider',
min=1,
max=25,
step=1,
value=1,
marks={i: f'{i} N' for i in range(1, 25, 1)},
included=False, # To remove the color fill
updatemode='drag', # Update slider continuously while dragging
tooltip={'placement': 'bottom'} # Show tooltip at the bottom
),
dbc.Row(
children=[
dbc.Col(
dcc.Graph(id='sampled-signal'),
width=6
),
dbc.Col(
dcc.Graph(id='sampled-fourier-signal'),
width=6
)
]
)
],
style={'backgroundColor': 'white', 'padding': '20px'}
)
# Callback to update graphs
@app4.callback(
[Output('sampled-signal', 'figure'),
Output('sampled-fourier-signal', 'figure')],
[Input('frequency-slider', 'value'),
Input('samples-slider', 'value'),
Input('N-slider','value')]
)
def update_graph(selected_frequencies,selected_samples,selected_N):
# Generate the signal with the selected frequencies
signal = np.zeros(len(t))
for f in range(selected_frequencies[0], selected_frequencies[1] + 1):
signal += np.sin(2 * np.pi * f * t)
# Correct downsampling approach
sample_rate = selected_samples # Directly using selected_samples as an integer
downsampling_factor = int(fs / sample_rate) # Downsampling factor
# Downsample the signal by selecting every nth sample
downsampled_signal = signal[::downsampling_factor]
# Recalculate L for the downsampled signal if needed
L_downsampled = len(downsampled_signal)
# Calculate the DFT of the downsampled signal
sampled_signal_fft = fft(downsampled_signal)
# Calculate magnitude of Fourier coefficients (amplitude) for the downsampled signal
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
n = len(downsampled_signal)
T = 1 / sample_rate # Recalculate the sampling period with the selected sample rate
# After calculating the magnitude of Fourier coefficients
sampled_magnitude = np.abs(sampled_signal_fft)[:L_downsampled // 2]
# Create a new frequency vector that includes the replicated frequencies
# First, calculate the original frequency bins for the positive frequencies
original_sampled_freqs = fftfreq(L_downsampled, T)[:L_downsampled // 2]
shifted_freqs_negative = original_sampled_freqs - 1/T # Shift for -1/Ts
shifted_freqs_positive = original_sampled_freqs + 1/T # Shift for +1/Ts
# Concatenate the original and shifted (replicated) frequencies and magnitudes
# This includes the original frequencies, and the replications at -1/Ts and +1/Ts
replicated_freqs = np.concatenate([shifted_freqs_negative, original_sampled_freqs, shifted_freqs_positive])
replicated_magnitude = np.concatenate([sampled_magnitude, sampled_magnitude, sampled_magnitude])
# Sort the replicated frequencies and magnitudes in ascending order for proper plotting
sorted_indices = np.argsort(replicated_freqs)
sorted_replicated_freqs = replicated_freqs[sorted_indices]
sorted_replicated_magnitude = replicated_magnitude[sorted_indices]
# Generate an array of sample indices for the downsampled signal
sample_indices = np.arange(len(downsampled_signal))
# Create data for the markers on top of the stems
sampled_signal_markers = go.Scatter(
x=sample_indices,
y=downsampled_signal,
mode='markers',
marker=dict(color='#00CC96', size=5),
name='Markers'
)
# Create a list to hold all the scatter plots for the stems
sampled_signal_stems = []
# For each sample point, create a line from the x-axis to the point
for x, y in zip(sample_indices, downsampled_signal):
sampled_signal_stems.append(
go.Scatter(
x=[x, x],
y=[0, y],
mode='lines',
line=dict(color='#00CC96', width=0.3),
showlegend=False # Hide legend for each stem
)
)
# Combine markers and stems for the plot
sampled_signal_data = [sampled_signal_markers] + sampled_signal_stems
# Update the figure with the new layout
sampled_signal_figure = {
'data': sampled_signal_data,
'layout': go.Layout(
title='Sampled Signal',
xaxis={'title': 'Sample Number', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Normalize the frequency values
normalized_freqs = original_sampled_freqs / (1/T)
# Since you previously concatenated and sorted for replication, ensure to apply normalization there as well
normalized_replicated_freqs = np.concatenate([
(shifted_freqs_negative / (1/T)), # Normalize -1/Ts shifted frequencies
normalized_freqs, # Already normalized original frequencies
(shifted_freqs_positive / (1/T)) # Normalize +1/Ts shifted frequencies
])
# Sort the normalized and replicated frequencies for proper plotting
sorted_indices = np.argsort(normalized_replicated_freqs)
sorted_normalized_replicated_freqs = normalized_replicated_freqs[sorted_indices]
sorted_replicated_magnitude = replicated_magnitude[sorted_indices]
# Let's assume N is defined; for example, N could be 10 for sampling every 1/10th
N=selected_N
# Adjust the sampling interval for selecting frequencies and magnitudes
sampled_indices = np.arange(0, len(sorted_normalized_replicated_freqs), N)
sampled_normalized_replicated_freqs = sorted_normalized_replicated_freqs[sampled_indices]
sampled_replicated_magnitude = sorted_replicated_magnitude[sampled_indices]
markers_data = go.Scatter(
x=sampled_normalized_replicated_freqs,
y=sampled_replicated_magnitude,
mode='markers',
marker=dict(color='#00CC96', size=5),
name='Markers',
showlegend=False
)
# Create a list to hold all the scatter plots for the stems
stems_data = []
# For each point, create a line from the x-axis to the point
for x, y in zip(sampled_normalized_replicated_freqs, sampled_replicated_magnitude):
stems_data.append(
go.Scatter(
x=[x, x],
y=[0, y],
mode='lines',
line=dict(color='#00CC96', width=0.25),
showlegend=False # Hide legend for each stem
)
)
# Combine markers and stems for the plot
data = [markers_data] + stems_data
# Update the figure with the new layout
sampled_fourier_transform_figure = {
'data': data,
'layout': go.Layout(
title='Sampled Normalized Fourier Transformation of Sampled Signal',
xaxis={'title': 'Sampled Normalized Frequency (f/fs)', 'showgrid': True, 'gridcolor': 'LightGrey'},
yaxis={'title': 'Amplitude', 'showgrid': True, 'gridcolor': 'LightGrey'},
margin={'l': 40, 'b': 40, 't': 40, 'r': 40},
hovermode='closest',
paper_bgcolor='white',
plot_bgcolor='white',
template='plotly_white'
)
}
# Return this additional figure in the callback
return sampled_signal_figure, sampled_fourier_transform_figure
if __name__ == '__main__':
app4.run_server(debug=True, port=8051)
Φασματική Ανάλυση#
Για τον υπολογισμό της ενέργειας ή ισχύος της κυματομορφής \(x(t)\), ανάλογα με την περίπτωση σήματος, ισχύει
όπου για σήματα ισχύος \(S_Χ(f)\) είναι η πυκνότητα φάσματος ισχύος (Power Spectral Density – PSD) της \(x(t)\). Για σήματα διακριτού χρόνου που προκύπτουν από δειγματοληψία της \(x(t)\) με περίοδο \(T_s\), οι αντίστοιχες σχέσεις υπολογισμό της ενέργειας ή ισχύος γίνονται
Ένας απλός τρόπος να εκτιμηθεί η πυκνότητα φάσματος ισχύος της κυματομορφής \(x(t)\) είναι να ληφθεί ο DTFT των δειγμάτων του σήματος και μετά να υψωθεί στο τετράγωνο το μέτρο του αποτελέσματος. Αυτός ο εκτιμητής αποκαλείται περιοδόγραμμα (periodogram). Το περιοδόγραμμα ενός πεπερασμένου μήκους \(L\) σήματος \(x[n]\) ορίζεται ως
όπου \(X_d(\phi)\) o DTFT του σήματος. Με το μήκος \(L\) να τείνει στο άπειρο, το περιοδόγραμμα \(P_{xx}(f)\) τείνει στην πυκνότητα φάσματος ισχύος \(S_Χ(f)\). Ο υπολογισμός του περιοδογράμματος σε πεπερασμένο πλήθος συχνοτήτων \(kf_s/Ν\), \(k=0\), \(1\), \(…\) , \(Ν\) δίνει
όπου \(X_k\) και ο DFT της πεπερασμένου μήκους \(L\) σειράς δειγμάτων του σήματος. Η ισχύς του σήματος είναι τότε
όπου η τελευταία ισότητα προκύπτει από το θεώρημα Parseval, που για την περίπτωση του DFT εκφράζεται ως:
Στην ειδική περίπτωση περιοδικών σημάτων έχουμε
όπου \(X[k]\) οι συντελεστές του αναπτύγματος σε σειρά Fourier και \(T_o\) η περίοδος του σήματος.
Μέρος 3: Εφαρμογή Α#
# Part 1: Create the signal
Fs = 1000 # Sampling frequency 1000 Hz
Ts = 1 / Fs # Sampling period
L = 1000 # Length of signal (number of samples)
T = L * Ts # Duration of signal
t = np.arange(0, (L - 1) * Ts, Ts) # Time vector
x = np.sin(2 * np.pi * 30 * t) \
+ 0.8 * np.sin(2 * np.pi * 80 * (t - 2)) \
+ np.sin(2 * np.pi * 60 * t) # 60 Hz component
# Time domain plot of x
fig = go.Figure()
fig.add_trace(go.Scatter(x=t, y=x, mode='lines', line=dict(color='#00CC96')))
fig.update_layout(title='Time domain plot of x', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='t (sec)',
yaxis_title='Amplitude',
template='plotly_white',
xaxis=dict(range=[0, 0.2]),
yaxis=dict(range=[-2.5, 2.5]))
fig.show()
# Fourier transform
def nextpow2(i):
n = 1
while n < i:
n *= 2
return n
N = nextpow2(L) # Length of Fourier transform
Fo = Fs / N # Frequency resolution
f = np.arange(0, N) * Fo # Frequency vector
X = np.fft.fft(x, N) # Compute DFT for N points
# Frequency domain plot of x (positive frequencies)
fig = go.Figure()
fig.add_trace(go.Scatter(x=f[1:N], y=abs(X[1:N]), mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Frequency domain plot of x', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='f (Hz)',
yaxis_title='Amplitude',
template='plotly_white')
fig.show()
# Shift frequencies to center
f = f - Fs / 2
X = np.fft.fftshift(X)
# Two-sided spectrum of x
f_shifted = f - Fs/2
fig = go.Figure()
fig.add_trace(go.Scatter(x=f_shifted, y=abs(X), mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Two sided spectrum of x', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='f (Hz)',
yaxis_title='Amplitude',
template='plotly_white')
fig.show()
# Calculate power
power = np.multiply(X, np.conj(X)) / N / L
# Periodogram
fig = go.Figure()
fig.add_trace(go.Scatter(x=f_shifted, y=power.real, mode='lines', line=dict(color='#1F77B4'))) # Use .real here if necessary
fig.update_layout(title='Periodogram', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='Frequency (Hz)',
yaxis_title='Power',
template='plotly_white')
fig.show()
# Part 2 Προσθέστε θόρυβο στο σήμα
# Συμπληρώστε τον κώδικα για τη δημιουργία του σήματος θορύβου n με τη βοήθεια της συνάρτησης randn.
# Το διάνυσμα θορύβου n θα πρέπει να είναι του ίδιου μεγέθους με αυτό της ημιτονοειδούς κυματομορφής x του πρώτου μέρους.
# Σχεδιάστε το σήμα θορύβου στο διάστημα από 0 έως 0.2 sec και κλίμακα σε από -2 έως 2.
# Υπολογίστε το περιοδόγραμμα του n και σχεδιάστε την πυκνότητα φάσματος ισχύος του σήματος θορύβου.
# Προσθέστε το σήμα θορύβου και το x για να λάβετε το σήμα με θόρυβο s.
# Σχεδιάσατε το σήμα με θόρυβο s στο πεδίο του χρόνου στην περιοχή 0 έως 0.2 sec
# και κλίμακα από -2 έως 2 καθώς και το αμφίπλευρο φάσμα του.
rand_n = np.random.randn(np.size(x))
fig = go.Figure()
fig.add_trace(go.Scatter(x=t, y=rand_n, mode='lines', line=dict(color='#00CC96')))
fig.update_layout(title='Time domain plot of n', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='t (sec)',
xaxis_range=[0, 0.2],
yaxis_title='Amplitude',
yaxis_range=[-2, 2],
template='plotly_white')
fig.show()
# Correction for N calculation using bitwise operator
N = 2^nextpow2(L)
Fo = Fs/N
f = (np.arange(0, N)) * Fo
f_shifted = f - Fs/2
rand_N = np.fft.fft(rand_n, N)
rand_N = np.fft.fftshift(rand_N)
power_n = np.multiply(rand_N, np.conj(rand_N)) / N / L
fig = go.Figure()
fig.add_trace(go.Scatter(x=f_shifted, y=power_n.real, mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Periodogram of n', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='Frequency (Hz)',
yaxis_title='Power',
template='plotly_white')
fig.show()
s = x + rand_n
fig = go.Figure()
fig.add_trace(go.Scatter(x=t, y=s, mode='lines', line=dict(color='#00CC96')))
fig.update_layout(title='Time domain plot of s', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='t (sec)',
xaxis_range=[0, 0.2],
yaxis_title='Amplitude',
yaxis_range=[-2, 2],
template='plotly_white')
fig.show()
S = np.fft.fft(s, N)
S = np.fft.fftshift(S)
fig = go.Figure()
fig.add_trace(go.Scatter(x=f_shifted, y=np.abs(S), mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Two sided spectrum of s', title_x=0.5, title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='f (Hz)',
yaxis_title='Amplitude',
template='plotly_white')
fig.show()
# Part 3. Πολλαπλασιασμός σημάτων
# Συμπληρώστε τον κώδικα δημιουργίας ενός ημιτονοειδούς σήματος συχνότητας
# 100 Hz και πολλαπλασιάστε με το προηγούμενο σήμα s.
# Τα δύο σήματα θα πρέπει να είναι του ίδιου μεγέθους.
# Σχεδιάστε το αποτέλεσμα στο πεδίο του χρόνου στην περιοχή 0 έως 0.2 sec
# και κλίμακα από -2 έως 2 καθώς και στο πεδίο της συχνότητας
# χρησιμοποιώντας τη συνάρτηση fftshift.
z=np.sin(2*np.pi*100*t)
y= np.multiply(z,s)
# Time domain plot of y
fig = go.Figure()
fig.add_trace(go.Scatter(x=t, y=y, mode='lines', line=dict(color='#00CC96')))
fig.update_layout(title='Time domain plot of y', title_x=0.5,
title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='t (sec)', xaxis_range=[0, 0.2],
yaxis_title='Amplitude', yaxis_range=[-2, 2],
template='plotly_white')
fig.show()
N = 2^nextpow2(L)
Fo = Fs/N
f = (np.arange(0, N)) * Fo - Fs/2
Y = np.fft.fft(y, N)
Y = np.fft.fftshift(Y)
# Frequency domain plot of y
fig = go.Figure()
fig.add_trace(go.Scatter(x=f, y=np.abs(Y), mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Frequency domain plot of y', title_x=0.5,
title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='f (Hz)', yaxis_title='Amplitude',
template='plotly_white')
fig.show()
Μέρος 4: Εφαρμογή Β#
Να γραφεί σε Python συνάρτηση φασματικής ανάλυσης, παρόμοια με την signal.welch(): θα δέχεται ως είσοδο διάνυσμα πραγματικού σήματος καθώς και τη συχνότητα δειγματοληψίας, \(F_s\), και θα σχεδιάζει τη μονόπλευρη φασματική πυκνότητα του σήματος στην περιοχή \([0-F_s/2)\). Το σήμα θα τεμαχίζεται σε τμήματα μήκους ίσου με τη δύναμη του \(2\) την πλησιέστερη στο \(1/8\) του συνολικού του μήκους, αλλά όχι μικρότερου από 256. Τα τμήματα θα είναι επικαλυπτόμενα κατά \(50\%\). Το τελευταίο τμήμα, εάν υπολείπεται σε μήκος των άλλων, θα αγνοείται. Θα υπολογίζεται με FFT το φάσμα κάθε τμήματος και θα λαμβάνεται η μέση τιμή όλων των τμημάτων. Η συνάρτηση να δοκιμαστεί με το σήμα του παραδείγματος 1.1 και να συγκριθεί το αποτέλεσμα με το αντίστοιχο της signal.welch().
def pwelch(x,Fs):
Ts=1/Fs
L=np.size(x)+1
T=L*Ts
N = 2^nextpow2(L)
Fo=Fs/N
f=np.arange(0,N)*Fo
window_size = nextpow2(np.size(x)/8)
if (window_size<256):
window_size=256
windows = np.size(x)//(window_size//2)-1
indexer = np.arange(window_size)[None, :] + (window_size//2)*np.arange(windows)[:, None]
windowed_x = x[indexer]
avg_pwr=0
for window in windowed_x:
window = window * np.hanning(np.size(window))
L=np.size(window)+1
T=L*Ts
N = 2^nextpow2(L)
Fo=Fs/N
f=np.arange(0,N)*Fo
window_fft=np.fft.fft(window,N)
power=np.multiply(window_fft,np.conj(window_fft))/N/L
avg_pwr=avg_pwr+power
avg_pwr=avg_pwr/windows
fig = go.Figure()
fig.add_trace(go.Scatter(x=f[:N//2], y=avg_pwr[:N//2].real, mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Periodogram pwelch()', title_x=0.5,
title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='Frequency (Hz)', yaxis_title='Power',
template='plotly_white')
fig.show()
return f[np.arange(0,N//2)], avg_pwr[np.arange(0,N//2)]
Fs=500
f1,Pxx1 = pwelch(x,Fs)
f2,Pxx2 = signal.welch(x,fs=Fs)
fig = go.Figure()
fig.add_trace(go.Scatter(x=f2, y=Pxx2, mode='lines', line=dict(color='#1F77B4')))
fig.update_layout(title='Periodogram signal.welch()', title_x=0.5,
title_font=dict(size=20, color='black', family="Arial, sans-serif"),
xaxis_title='Frequency (Hz)', yaxis_title='Power',
template='plotly_white')
fig.show()